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I.  INTRODUCTION 


The  performance  characteristics  of  a  hydrodynamic  propulsor  are 
determined,  to  a  large  extent,  by  complex  three-dimensional  viscous  flow  in  the 
blade  rows.  The  capability  to  compute  the  complex  flow  field  for  use  in  the 
design  and  analysis  of  performance  of  a  propulsor  would  be  a  valuable  tool  for 
improving  performance. 

Research  under  the  present  contract  addressed  computation  of  the  steady 
flow  field  in  shrouded  propulsor  blade  rows.  A  prime  consideration  in  the 
choice  of  a  computation  procedure  is  the  multiplicity  of  length  scales  in  the 
shrouded  propulsor  flow  field  that  need  to  be  adequately  resolved.  The 
differing  dominant  flow  mechanisms  at  the  hub  and  tip  of  the  blade,  blade 
boundary  layers,  and  a  core  inviscid  flow  characterize  the  multiple  length 
scales.  A  single  computation  procedure,  useful  in  a  design  and  analysis 
context,  that  adequately  resolved  all  of  the  flow  length  scales  was  deemed 
uneconomical  due  to  the  large  number  of  grid  points  that  would  be  necessary  for 
the  computation.  The  approach  pursued,  under  the  present  contract,  was  to 
study  the  flow  field  in  the  tip  and  hub  regions  of  the  blade  in  separate  more 
tractable  computations  (zonal  approach). 

The  dominant  flow  mechanisms  in  the  blade  tip  region  of  a  shrouded 
propulsor  are  due  to  tip  clearance  effects  and  the  shroud  boundary  layer.  A 
preliminary  analysis  of  the  physics  of  the  shrouded  propulsor  tip  flow  field 
indicated  that  a  unique  spatial-marching  procedure  developed  at  SRA  could  be 
utilized  to  compute  this  flow  field.  The  primary  motivation  in  adopting  this 
approach  was  the  very  high  density  grids  required  for  adequate  resolution  of 
the  complex  flow  structure  of  the  shrouded  tip  flow  field.  Spatial-marching 
procedures  can  utilize  these  high  density  grids  far  more  economically  than 
methods  for  the  three-dimensional  Navier-Stokes  equations.  The  accuracy  and 
high  computational  efficiency  of  the  spatial-marching  analysis  has  been 
demonstrated  in  applications  to  internal  and  external  flows  [1-4].  The  run 
time  needed  for  solution  is  10  to  100  times  less  than  that  required  for  the 
Navier-Stokes  equations.  Success  of  the  proposed  approach  could  provide  the 
capability  to  address  a  very  difficult  viscous  flow  problem  with  both  economy 
and  accuracy. 

The  approximate  flow  equations  solved  in  the  spatial-marching  procedure 
were  judged  to  be  inadequate  for  the  hub  region  of  the  propulsor  flow  field. 


"  N 

/ 

Uist 

/)'/ 

-1- 


This  is  due  to  the  fact  that  the  approximations  made  are  inaccurate  in  the 
region  of  the  blunt  leading  edge  of  the  propulsor  blade  at  the  hub. 

Furthermore,  the  hub  flow  field  often  contains  a  large  region  of  separated  flow 
at  the  corner  of  the  suction  surface  of  the  blade  and  the  hub.  Large  regions 
of  streamwise  flow  separation  are  beyond  the  scope  of  the  spatial-marching 
procedure.  Solution  of  the  Navier-Stokes  equations  was  proposed  for  the 
propulsor  hub  flow  field. 

The  present  report  describes  progress  made  under  the  Propulsor  Blade  Flows 
contract  towards  computations  of  the  tip  and  hub  flow  fields  of  a  shrouded 
propulsor.  The  next  section  of  the  report  describes  the  application  of  the 
spatial-marching  procedure  to  the  shrouded  propulsor  blade  tip  flow  field.  The 
third  section  of  the  report  describes  progress  towards  application  of  a 
Navier-Stokes  solution  procedure  to  the  hub  flow  field. 

2.  SHROUDED  PROPULSOR  BLADE  TIP  FLOW  FIELD 

An  analysis  of  the  application  of  the  spatial-marching  procedure  to  the 
shrouded  tip  flow  field  was  carried  out  for  the  purpose  of  defining 
approximations  and  boundary  conditions.  Such  an  analysis  also  provided  insight 
into  the  physical  nature  of  the  flow  field  in  the  region  of  the  shrouded  tip. 

A  model  problem  was  devised  for  this  purpose  and  consisted  of  a  constant 
thickness  blade  with  a  square  tip  and  a  moving  shroud  wall  (Figure  A).  The 
moving  shroud  wall  simulates  relative  motion  between  the  propulsor  blade  tip 
and  the  shroud.  This  simplified  problem  geometry  was  chosen  to  focus  efforts 
on  key  elements  of  the  shrouded  tip  flow  structure  affected  by  the  physical 
approximations  made  in  the  primary/secondary  flow  equations  solved  by 
spatial-marching. 

A  critical  aspect  of  computing  the  shrouded  tip  flow  field  and  the  flow 
field  of  the  model  problem  is  the  mass  flux  through  the  tip  clearance  region. 
This  mass  flux  is  controlled,  to  a  large  extent,  by  transverse  and  streamwise 
pressure  gradients  in  the  tip  clearance  region,  viscous  effects  of  the  relative 
motion  between  the  shroud  and  the  blade  tip,  and  the  mass  flux  through  the 
lateral  boundaries  of  the  computational  domain.  Approximating  streamwise 
pressure  gradients  as  a  means  to  obtain  equations  that  are  well-posed  for 
spatial-marching  was  judged  to  be  inappropriate  as  a  physical  approximation  due 
to  the  strong  influence  of  this  pressure  gradient  on  the  computed  mass  flux 
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through  the  tip  clearance.  Typically,  approximate  streamwise  pressure 
gradients  are  obtained  from  the  potential  flow  through  the  geometry  under 
consideration.  For  the  model  problem  the  potential  flow  streamwise  pressure 
gradient  is  zero.  Computations  for  the  model  problem  with  the  streamwise 
pressure  gradient  set  to  zero  produced  non-physical  flow  fields. 

Generalized  primary/secondary  flow  equations  have  been  derived  at  SRA 
(1,  2],  well-posed  for  a  spatial-marching  procedure,  that  introduce  no 
approximations  for  pressure  gradients.  Approximations  are  introduced,  instead, 
on  the  scalar  potential  component  of  a  decomposed  transverse  velocity  field. 

The  scalar  potential  velocity  component  is  defined  to  be  normal  to  a  local 
primary  flow  direction.  Typically,  the  primary  flow  direction  is  chosen  to  be 
the  direction  of  the  potential  flow  streamline  at  the  point.  Appendix  1 
describes  the  generalized  pr imary/secondary  flow  equations  and  a  sequentially 
decoupled  implicit  algorithm  used  to  solve  the  equations. 

Application  of  the  generalized  primary/secondary  flow  equations  to  the 
shrouded  tip  flow  field  required  integrating  the  spatial-marching  procedure 
with  the  treatment  of  the  blade  tip  region  as  a  local  zone.  The  model  problem 
(Figure  A)  provided  a  simple  geometry  to  test  developed  procedures.  A  key 
element  in  this  development  was  to  construct  appropriate  boundary  conditions  at 
lateral  boundaries  for  the  generalized  primary/secondary  flow  equations. 
Accuracy  of  the  approach  could  be  judged  by  the  magnitude  of  the  computed 
scalar  potential  velocities  in  the  flow  field.  Boundary  conditions  along  the 
other  boundaries  of  the  model  problem  have  been  developed  in  previous  work  or 
were  straightforward.  Boundary  conditions  on  the  shroud  and  blade  surfaces 
were  derived  from  the  no-slip  condition  on  a  solid  surface.  Inboard  boundary 
conditions  have  been  developed  for  the  tip  vortex  flow  field  [3]  and  are 
suitable  for  the  shrouded  tip  flow  field.  However,  for  the  preliminary 
computations  carried  out,  the  flow  field  was  assumed  to  be  symmetric  along  the 
inboard  boundary.  These  boundary  conditions  simplified  the  development  and 
testing  of  lateral  boundary  conditions  appropriate  for  the  zonal  approach. 

Two  approaches  were  attempted  in  the  development  of  lateral  boundary 
conditions  for  the  primary/secondary  flow  equations  along  with  the  zonal 
approach.  These  conditions  needed  to  be  specified  in  terras  of  the  scalar  and 
vector  potentials  for  the  decomposed  transverse  velocity  field.  The  critical 
element  in  these  conditions  is  the  value  of  the  vector  potential  on  the  shroud 
that  represents  the  mass  flux  entrained  by  the  relative  motion  between  the 
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blade  tip  and  the  shroud.  The  first  approach  determined  the  value  of  the 
vector  potential  on  the  shroud  from  a  local  boundary  layer  analysis  at  the 
lateral  boundary.  However,  computations  for  the  model  problem  showed  large 
scalar  potential  velocities  and,  therefore,  inaccurate  solutions.  Small  errors 
in  mass  flux,  due  to  the  approximate  boundary  layer  analysis,  at  the  lateral 
boundary  were  magnified  at  the  small  tip  clearance  gap.  The  second  approach 
attempted  a  local  analysis  in  the  clearance  gap,  based  on  small  normal 
velocites  in  the  gap,  to  determine  the  vector  potential  on  the  shroud.  The 
second  approach  was  also  unsuccessful  due  to  inaccurate  solutions  from  a  large 
scalar  potential  velocity  field. 

All  of  the  lateral  boundary  conditions  developed  and  tested  in  the  present 
investigation  led  to  large  scalar  potential  velocities  near  the  tip  clearance 
gap.  Consequently,  the  physical  approximations  required  to  effect  a 
spatial-marching  solution  were  inaccurate  for  this  flow  application.  From 
these  results  the  conclusion  was  reached  that  considerable  further  work  and 
analysis  for  lateral  boundary  conditions  would  be  required  before  attempting  a 
zonal  approach  together  with  the  spatial-marching  procedure  for  computing  the 
shrouded  tip  flow  field. 

3.  A  COORDINATE  SYSTEM  FOR  THE  NAVIER-STOKES 
CALCULATION  OF  A  PROPELLER  HUB  FLOW  FIELD 

A  coordinate  sytem  is  one  component  in  the  application  of  finite- 
difference  numerical  techniques  to  fluid  dynamics  problems.  A  good  coordinate 
system  should  minimize  the  error  in  the  governing  equations  for  the  problem 
being  analyzed.  Failure  to  adequately  resolve  changes  in  the  fluid  dynamic 
properties  between  mesh  points  is  one  source  of  error.  Excessive  changes  in 
grid  spacing  between  neighboring  grid  points  as  well  as  excessive  coordinate 
curvature  and  non-orthogonality  also  produce  error  terms.  An  accurate 
calculation  requires  a  coordinate  system  with  adequate  resolution  of  the 
gradients  of  the  flow  variables  and  smoothness  of  the  computational  grid. 

In  generating  the  propeller  hub  coordinate  system,  care  has  been  taken  to 
resolve  near  wall  regions  and  to  provide  smoothness  of  the  grid  lines  and  of 
the  grid  spacing  along  each  line.  This  is  accomplished  by  building  a 
three-dimensional  coordinate  system  from  a  series  of  well-behaved  two- 
dimensional  coordinate  systems.  The  two-dimensional  coordinates  are 
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generated  to  have  favorable  smoothness  properties  and  the  subsequent  mappings 
are  designed  to  preserve  the  favorable  smoothness  properties  of  the 
coordinates . 

An  H-type  coordinate  system  about  an  uncambered  Joukowski  airfoil  is 
generated  using  potential  lines  and  streamlines  about  the  airfoil  at  zero 
incidence.  The  streamlines  are  mapped  to  horizontal  lines  at  a  specific 
distance  from  the  airfoil.  This  coordinate  system  is  shown  in  Figure  la.  A 
series  of  these  airfoils  can  be  stacked  to  form  a  cascade  of  zero  pitch.  This 
coordinate  system  provides  nearly  orthogonal  coordinates  and  fine  resolution 
near  the  blades,  as  shown  in  Figure  lb.  Grid  clustering  is  provided  in  the 
streamwise  direction  near  the  leading  and  trailing  edges.  Stacking  this 
coordinate  system  provides  continuity  of  location,  slope  and  curvature  for  the 
vertical  coordinate  lines  across  the  upper  and  lower  boundaries. 

Pitch  is  added  to  the  cascade  by  performing  a  shearing  transformation 
shown  in  Figure  2.  The  shearing  transformation  produces  a  midline  shape  of  a 
circular  arc  with  upstream  and  downstream  straight  line  extensions.  The 
circular  arc  need  not  coincide  with  the  blade  chord.  The  inflow  and  outflow 
boundaries  are  straight  vertical  lines.  The  shearing  transf ormation  distorted 
the  blade  shapes  and,  since  Joukowski  thickness  distributions  are  not  used  in 
practice,  the  desired  blade  thickness  distribution  is  imposed.  A  blade 
thickness  distribution  is  imposed  normal  to  the  raidline  of  the  sheared  blade. 
The  coordinate  lines  emanating  from  the  blade  are  relaxed  to  the  new  blade 
surface.  The  results  of  this  process  is  shown  in  Figure  2  for  a  blade 
thickness  distribution 
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where  X  is  the  distance  along  the  chord,  Y  is  the  distance  from  the  chord  to 
the  blade  and  C  is  the  chord  length. 

This  coordinate  system  can  be  stacked  to  form  a  cascade  of  airfoils  as 
shown  in  Figure  3.  At  the  upper  and  lower  boundaries  the  coordinates  are 
continuous  in  location,  slope  and  curvature  making  application  of  periodic 
boundary  conditions  very  natural.  In  addition,  the  coordinates  are  clustered 
vertically  near  the  blades  and  horizontally  near  the  leading  and  trailing  edges 
to  ensure  resolution  of  the  large  gradients  of  the  flow. 
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Specifying  blade  shape  and  coordinate  system  length  and  width,  coordinate 
systems  such  as  shown  in  Figure  2  are  mapped  into  9,z  surfaces  in  the  hub 
coordinate  system.  Note  that  in  these  surfaces  the  radius  will  be  a  function 
of  the  z  coordinate  as  shown  in  the  side  view  of  the  hub  in  Figure  4.  These 
r,z  coordinates  are  computed  as  the  leading  edge  of  a  thick  airfoil  (turned 
backwards)  using  potential  lines  and  streamlines  as  outlined  above. 

The  x,y  coordinates  of  Figure  2  are  mapped  into  the  r,z  surfaces  of  Figure 
4.  Figure  5  shows  three  blades  on  a  hub  with  the  r,z  coordinate  lines  of 
Figure  4  shown  in  white.  Figure  6  shows  the  inner  and  outer  8,z  surfaces. 

These  blades  are  stacked  on  their  leading  edges,  although  the  choice  of 
stacking  location  can  be  specified.  Blade  twist  is  seen  from  the  twist  of  the 
upstream  extension  of  the  coordinate  system.  Five  regions  of  fine  resolution 
are  contained  in  this  coordinate  system:  1)  hub  surface,  2)  hub  trailing  edge, 
3)  blade  surfaces,  4)  blade  leading  edge  and  5)  blade  trailing  edge. 
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CONCLUSIONS 


(a)  Application  of  a  spatial-marching  procedure  to  compute  the  shrouded  tip 
flow  field  was  unsuccessful  due  to  the  inability  to  construct  boundary 
conditions  compatible  with  the  approximations  in  the  flow  equations. 

(b)  A  grid  generation  procedure  for  the  propulsor  hub  flow  field  to  construct 
computational  grids  for  use  with  a  Navier-Stokes  solver  was  successfully 
completed. 
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Appendix  I 


1.  Background 


The  pr imary/secondary  flow  equations  derived  in  Ref.  1  are  an  approximate 
set  of  flow  equations  suitable  for  complex  three-dimensional  flow  problems  with 
a  dominant  flow  direction  that  could  contain  large  streamwise  vorticity. 
Streamwise  viscous  diffusion  is  small  in  such  flows,  compared  to  transverse 
diffusion,  and  is  neglected  in  the  primary/secondary  flow  equations.  The  small 
scalar  potential  approximation  makes  the  primary/secondary  flow  equations  a 
well-posed  initial-value  problem  in  the  streamwise  coordinate  direction.  The 
streamwise  coordinate  is  nominally  aligned  with  the  dominant  flow  direction.  An 
important  element  in  the  derivation  of  the  primary/secondary  flow  equations  is 
that  the  approximations  are  associated  with  directions  identified  by  a  vector 
basis  at  each  point  in  the  flow  field,  distinct  from  coordinate  directions.  The 
vector  basis  for  the  approximations  provides  considerable  flexibility  in 
choosing  convenient  coordinate  systems  for  complex  three-dimensional  geometries, 
while  at  the  same  time  making  valid  flow  approximations. 

The  motivation  in  deriving  the  primary/secondary  flow  equations  is  to 
obtain  a  set  of  approximate  equations  that  can  be  solved  economically  compared 
to  solving  the  full  Navier-Stokes  equations.  The  primary/secondary  flow 
equations  can  be  solved  as  an  initial-value  problem  in  the  streamwise  coordinate 
direction,  compared  to  the  boundary-value  Navier-Stokes  equations  in  the  same 
direction,  with  potential  savings  in  computational  effort  and  storage 
requirements.  These  savings  provide  the  ability  to  solve  a  large  number  of  flow 
cases,  with  sufficient  grid  points  to  resolve  complex  three-dimensional  flow 
fields,  from  limited  resources.  An  algorithm  is  developed  to  solve  the 
pr imary/secondary  equations  that  provides  an  order  of  magnitude  (factor  of  10  to 
20)  advantage  in  computational  effort  over  an  efficient  Navier-Stokes  solver. 

The  next  section  of  this  Appendix  contains  a  list  of  the  primary/secondary 
flow  equations  derived  in  Ref.  1.  Properties  of  the  equations  that  affect  the 
choice  of  the  solution  algorithm  are  discussed.  A  Sequentially  Decoupled 
Implicit  algorithm  is  developed  as  an  efficient  means  of  solving  the 
pr imary/secondary  equations.  Numerical  procedures  that  comprise  the  algorithm 
and  their  implementation  are  discussed. 
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2.  Primary/Secondary  Flow  Equations 

The  primary /secondary  flow  equations  are  listed  here  for  convenience  of 
presenting  the  solution  algorithm.  Derivation  of  these  equations  can  be  found 
in  Ref.  1.  In  the  equations  that  follow,  ik  denotes  the  unit  vector  basis 
associated  with  the  orthogonal  reference  line  coordinate  system.  The  vector  ij 
is  aligned  with  the  reference  line  direction  (x^)  along  which  the  primary/ 
secondary  flow  equations  are  an  initial-value  problem.  The  vector  basis  ek 
and  the  reciprocal  basis  e^  are  defined  for  the  purpose  of  approximating 
flow  equations. 

ck*e£  =  ekme,nf-  -  <5k  (2.1) 

Th-->  flow  field  velocity  vector  is  denoted  by  U, 

U  =  u1i1  +  U2 i 2  +  u  3 i 3  =  Upe 1  +  ve2  +  we3  (2.2) 

The  velocity  component  Up  is  associated  with  the  primary  flow  direction  and 
components  v  and  w  with  the  secondary  flow  field.  The  secondary  flow  velocity 
vector  V  can  be  defined  as 


V  =  v,2  +  we 3  (2.3) 

The  secondary  flow  velocity  components  are  decomposed  into  scalar  potential 
(v,.. ,  Wq)  and  vector  potential  components  (v^,  w^) 

V  =  +  v^ 

W  =  +  w^  (2.4) 

V  =  v*  +  v* 

The  scalar  potential  velocity  components  can  be  defined  in  terms  of  a  scalar 
potential  (<J>)  and  the  vector  potential  components  in  terms  of  a  vector 
oo tent  ial  (  >J»)  • 
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=  (e^V^ej  )i* 

£,  =  2,3  (no  sum)  (2.5) 


—  £  T  —  —  — 

VgX^i!  =  (pe  JV^-ej)i£ 


where  p  is  density  and  7S  is  the  surface  operator  defined  as 


7  =  7  -  V 

s  p 


with 


7  =  i i  (i-i  -  V  ) 

P 


(2.6) 


All  variables  in  the  flow  equations  are  assumed  to  be  nondimensional ,  having 
been  normalized  by  a  suitable  combination  of  reference  quantities.  The 
primary/secondary  flow  equations  are  written  for  the  dependent  variables, 
streamwise  velocity  Up,  scalar  potential  4>  >  vector  potential  i|>,  vorticity  ft , 
pressure  p,  total  enthalpy  E,  and  density  p. 


Scalar  potential  equation  (continuity): 


7  *(p7  i)  +  7  •  (pe^V*  e  .  )ii  +  V«  (pe^u  )i  =  0  (2.7) 

s  s  p  j  p  k 


Streamwise  momentum  equation: 


p  'J  •  7  u 

P 


pO- (0*7 )ej 


+  e^ •  Vp 


ej,7xyVxU  =  0 


(2.8) 


where  Re  is  the  Reynolds  number  defined  trora  reference  quantities,  u  the 
normalized  coefficient  ot  molecular  viscosity,  and  7  is  an  approximation  to  7 
in  which  second  derivatives  along  the  ij  coordinate  have  been  neglected. 


Vector  potential  equation  (definition  of  vorticity): 


V  xi.  (7  x^i.)  -  ftii  =  0 
s  p  s 


1 


(2.9) 


Vorticity  transport  equation: 


puD  3Q  -  fk  f.  -  pun  aP£k  -  - 

-f c  nr  *  's'0™  +  V(<!  VD  ut  +  (5 s-ir ^  <v^)i* 


+  WsPV)«[vs,tk(S-ik)]i1[  +  *1  Ve'‘vV"!1‘>  -  o 


(2.10) 
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where  C.  and  D  are  defined  in  (lj  and  V  =  gV  +  V, 

k  ijj 

values  of  0  or  1,  6  =  0  identifies  the  small  scalar 


.  The  parameter  g  takes  on 
potential  approximation. 


Pressure  equation: 


Ik-  -  tk  — 

Vs.(o  e^’VpH^  +  Vs-(e  Fk. )i ^  =  0  £  =  2,3  (no  sum)  (2.11) 


in 

F, 

< 


is  a  combination  of  terms  from  transverse  momentum  equations  and  is  defined 
[1].  The  secondary  flow  velocity  vector,  V,  is  written  as  V  =  gV^+  in 
and  6=0  identifies  the  small  scalar  potential  approximation. 


Energy  equation: 


oU*7E  -  (RePr)-1  7><V  E  -  (  —  1)  M2Re_1  T  u$  -  pU-V  xpV  xq 

s’  r  s  s 

(2.12) 

-  Pr_1  VgicV^  | q  | 2  / 2 ]  =  0 


where  q  =  u  e1  +  (Bv,  +  v,  )e2  +  (gw,  +  w,  )e3  , 

P  <t>  ^  4>  ’J' 

$  is  the  dissipation  function,  <  the  coefficient  of  thermal  conductivity,  Pr 
the  Prandtl  number,  and  Mr  the  reference  Mach  number.  Simpler  forms  of  the 
energy  equation  that  can  be  used  at  low  and  moderate  Mach  numbers  are  described 
in  [  1  ] . 

State  equation: 


p  =  pE/yM2  -  (y — 1 )p  q  2 /2y 


(2.13) 
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Equations  (2.7)  to  (2.13),  along  with  the  definition  of  the  velocity 
decomposition  (2.5),  are  the  complete  set  of  primary/secondary  flow  equations. 
The  small  scalar  potential  approximation  (8=0)  makes  these  equations  a 
well-posed  initial  value  problem  along  the  streamwise  coordinate  (xj). 

3.  Solution  of  Primary/Secondary  Flow  Equations 

The  choice  of  a  solution  algorithm  for  the  primary/secondary  flow 
equations  must  exploit  the  small  scalar  potential  approximation  for  economy. 

The  potential  for  economy  lies  in  the  initial-value  nature  of  the  primary/ 
secondary  flow  equations  and  the  ability  to  solve  such  a  problem  by  a  spatial 
marching  scheme.  Such  a  scheme  will  be  developed  in  this  section  for  the 
primary/secondary  flow  equations. 

An  evolution  equation  can  be  written  in  the  form 

A  +  D(  4>)  =  0  (3.1) 

9t 

where  )  is  a  vector  of  dependent  variables,  t  is  a  time-like  variable  (e.g. 
streamwise  coordinate),  A  is  a  matrix  of  coefficients,  and  D  is  a  multi¬ 
dimensional  differential  operator.  Efficient,  implicit,  noniterative  marching 
schemes  have  been  developed  to  solve  equations  of  the  type  (3.1)  when  the 
matrix  A  is  well-conditioned.  An  example  of  such  a  scheme  is  the  LBI  scheme 
[5]  for  nonlinear  equations  and  is  based  on  the  Douglas-Gunn  ADI  scheme. 
However,  the  matrix  A  for  the  primary/secondary  flow  equations  is  singular  due 
to  the  lack  of  streamwise  derivatives  in  some  of  the  governing  equations.  This 
is  also  true  when  the  same  equations  are  written  in  terras  of  primitive 
variables  [1].  Coupled  implicit  schemes  such  as  the  LBI  scheme  require 
iteration  when  A  is  singular.  This  penalty  in  economy,  due  to  iteration  of  a 
fully  coupled  system  at  each  step,  is  avoided  in  the  present  sequentially 
decoupled  algorithm. 

The  problem  with  the  singularity  of  the  coefficient  matrix  of  the 
streamwise  derivative  is  circumvented  by  decoupling  the  individual  equations  of 
the  system  without  a  streamwise  derivative  and  solving  them  separately  at  each 
streamwise  step  of  a  spatial  marching  scheme  for  the  remaining  equations. 
However,  any  severe  penalty  in  economy  due  to  loss  of  stability  in  the  overall 
scheme  from  the  decoupling  would  be  unacceptable.  The  present  sequentially 
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decoupled  implicit  algorithm  decouples  the  primary /secondary  flow  equations 
into  subsystems  through  sequencing  and  linearization  of  the  subsystems.  Each 
of  the  subsystems  can  then  be  solved  by  economical,  stable,  implicit 
procedures. 

3.1  Sequentially  Decoupled  Implicit  Algorithm 

The  solution  algorithm  for  the  primary/secondary  flow  equations  is  developed  by 
considering  a  more  general  case  of  a  partial  differential  equation  in  the  form 

Mq  =  0  (3.2) 

M  is  a  matrix  operator  (nonlinear)  whose  elements  are  the  coefficients  and 

differential  operators  of  the  individual  equations  of  the  system  and  q  is  a 

vector  of  dependent  variables.  Further,  the  elements  of  q  are  assumed  to  be 

ordered  and  the  sequence  of  equations  in  M  arranged  on  the  basis  of  an 

association  between  an  equation  and  a  dependent  variable.  This  association 
between  individual  equations  of  the  system  and  dependent  variables  will  be 
clarified  with  development  of  the  algorithm.  The  matrix  M  is  decomposed  and 
(3.2)  written  in  the  form 


(L  +  D  +  u)q  =  0  (3.3) 

where  L,  D,  U  are  matrix  operators  containing  the  lower  diagonal,  diagonal  and 
upper  diagonal  elements  of  M,  respectively.  The  matrix  operators  L,  D  and  U 
are  linearized  by  the  following  linearization  scheme 

Ln’n+1  :  (qn+1,  qn)  k  =  1,  .  .  .  (j-l)  (3.4) 


where  t ( j  are  the  elements  of  L  and  n  is  the  marching  coordinate  index.  The 
operators  D  and  U  can  be  linearized  in  the  same  manner  as  L.  Equation  (3.3) 
can  be  written  as  a  linear,  decoupled  set  of  equations 

(Dn’n+1  +  Ln’n+1)  qn+1  +  Un’n+1  q"  =  0  (3.3) 
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or 


^n,n+l  n+1  ,.n,n+l  n+1  ,  „n,n+l  n, 

D  q  =  -  ( L  q  +  U  q) 


(3.6) 


Equation  (3 is  a  decoupled  system  of  equations  that  can  be  solved  by 
inverting  the  linear  diagonal  operators  in  sequence  using  implicit  schemes. 
Clearly,  the  success  of  the  sequentially  decoupled  algorithm  (3.6)  depends  on 
the  system  of  equations  (3.2)  and  the  choice  of  the  diagonal  matrix  operator, 

D.  The  diagonal  operator,  D,  is  determined  by  the  association  made  between 
individual  equations  of  the  system  and  dependent  variables. 

The  association  of  equations  with  variables  is  straightforward  for  the 
pr imary/secondary  flow  equation;  ft  with  (2.10),  i|>  with  (2.9),  p  with  (2.11), 

Up  with  (2.8),  $  with  (2.7),  and  E  with  (2.12).  For  simplicity,  the  implicit 
operator  (Dn,n  +  Ln,n  )  in  (3.5)  for  the  incompressible  form  of  the  primary/ 
secondary  flow  equations  can  be  written  as 


(Dn’n+1  +  Ln’n+1)  qn+1 


-  - - 

—  — 

dL  i  *  * 

n 

1  V2 

<P 

s 

*31  *32  V2  * 

p 

„  „  1  3 

0  0  Fk 

u 

p 

o  o  o  _!  v2 

3x  s 

♦ 

_  _ 

n+1 


(3.7) 


where  dllt  d,^  =  U*V  -  Re 


The  in  (3.7)  is  used  to  indicate  non-zero 
,n,n+l 


elements  of  the  explicit  operator,  U  *  .  The  vector  basis  is  assumed  aligned 

with  coordinate  directions  in  (3.7).  The  vector  basis  introduces  more  complex 
differential  operators  in  individual  equations  but  the  implicit  operator 
retains  essentially  the  same  form  as  (3.7).  Although  no  stability  analysis  of 
(3.6)  has  been  performed,  stability  of  (3.6)  for  the  primary/secondary  flow 
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equations  has  been  established  through  extensive  numerical  computation  of 
numerous  test  cases.  Further,  although  the  vorticity  and  vector  potential 
equations  (first  two  equations  in  (3.7))  can  be  uncoupled  in  (3.7),  they  are 
strongly  coupled  at  a  no-slip  boundary  by  boundary  conditions.  Consequently, 
the  vorticity  and  vector  potential  equations  are  solved  as  a  coupled  subsystem 
of  equations  to  avoid  stability  problems.  The  following  sequentially  decoupled 
subsystem  equations  are  obtained  for  the  compressible  primary/secondary  flow 
equations 

1)  vorticity  (2.10),  vector  potential  (2.9)  for  ft,  <J> 

2)  pressure  equation  (2.11)  for  p 

3)  energy  equation  (2.12)  for  E 

4)  streamwise  momentum  (2.8)  for  Up 

5)  equation  of  state  (2.  13)  for  p  (algebraic) 

6)  scalar  potential  equation  (2.7)  for  cf> 

The  sequence  completes  the  solution  of  the  primary/secondary  flow  equations  for 
one  streamwise  step  (n+1)  from  the  solution  at  the  previously  computed  step  (n) 
(initial  conditions  for  n=l).  Each  of  the  individual  subsystem  equations  can 
be  solved  by  a  convenient,  efficient  procedure.  Numerical  procedures  used  in 
the  present  work  for  solving  decoupled  subsystem  equations  are  described  in  a 
later  section. 

3.2  Differencing  Procedures  and  Boundary  Conditions 

The  governing  equations  (2. 7)— (2. 12)  are  transformed  to  a  body-fitted 
coordinate  system  in  which  the  equations  are  discretized.  Body-fitted 
coordinate  transformations  suitable  for  the  primary/secondary  flow  equations 
are  discussed  in  [1).  The  discretization  employs  two-point  backward 
differences  for  streamwise  derivatives  and  three-point  central  differences  for 
transverse  derivatives.  The  body-fitted  coordinate  transformation  allows  a 
variable  streamwise  step  size  and  concentration  of  grid  points  in  regions  of 
expected  large  gradients.  A  centered  Crank-Nico Ison  formulation  could  be  used 
for  streamwise  derivatives  but  the  formal  accuracy  of  the  scheme  would  not  be 
improved  because  of  the  first-order  accurate  linearization  scheme  (3.4)  used  to 
decouple  subsystem  equations. 
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Implicit  boundary  conditions  are  used  in  all  subsystem  equations. 

Boundary  conditions  that  are  not  straightforward  are  no-slip  conditions  on  a 
solid  boundary  and  boundary  conditions  on  pressure.  These  conditions  are 
described  here.  The  streamwise  velocity  Up  is  specified  as  zero  in  the 
streamwise  momentum  equation  (2.8)  on  a  no-slip  boundary.  The  condition 
v  =  w  =  0  must  be  expressed  in  terms  of  <J> ,  ifi  and  £1.  The  normal  component  of 
the  scalar  potential  velocity  is  specified  as  zero  by  in»V^  =  0,  where  i^  is 

the  unit  normal  vector  to  the  boundary.  V  is  specified  in  terms  of  <{>  to 

v 

obtain  a  Neumann  condition  for  the  scalar  potential  equation  (2.7).  The  normal 
component  of  the  vector  potential  velocity  is  specified  zero  by  prescribing  = 
0  on  the  boundary.  The  tangential  component  of  the  transverse  velocity  is 
specified  at  the  no-slip  boundary  by 


vl(v6  +  Ve  +  (wd, +  Ve  1 


3 1  - 


0 


(3.8) 


where  ic  is  the  unit  tangent  vector  to  the  boundary.  Equation  (3.8)  combined 
with  equation  (2.9)  is  used  to  obtain  £1  on  the  boundary  in  terms  of  <j>  and  \J). 
This  relation  for  £2,  along  with  prescribing  <J>  =  0,  provides  coupled  boundary 
conditions  for  the  vector  potential  (2.9)  and  vorticity  equations  (2.10).  The 
coupled  boundary  conditions  require  the  equations  to  be  solved  as  a  coupled 
subsystem  from  stability  considerations. 

Boundary  conditions  for  pressure  are  obtained  by  computing  the  normal 
gradient  of  pressure  on  the  boundary  from  transverse  momentum  equations.  The 
prescription  of  the  normal  gradient  of  pressure  in  the  pressure  equation  (2.11) 
results  in  a  Neumann  problem.  Solution  of  the  Neumann  problem  leaves  an 
undetermined  constant  in  the  transverse  plane  that  could,  however,  be  a 
function  of  the  streamwise  coordinate.  Pressure,  p,  is  written  as 

p  =  p  (xi )  +  p  (x^,  x2 ,  x3)  (3.9) 

03  V 


where  pm  is  the  undetermined  constant  of  the  Neumann  problem  and  pv  is  the 

solution  of  the  pressure  equation.  The  value  of  pm  is  required  in  the 

solution  of  the  streamwise  momentum  equation  (2.8)  in  the  form  of  the 

streamwise  gradient,  .  In  external  flows,  pm  can  be  determined  by  equating 

3  x  i 
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the  right-hand  side  of  (3.9)  to  prescribed  pressure  at  a  point  in  freestream  in 
each  transverse  plane.  Such  a  condition  is  not  available  in  internal  flows  and 
pm  has  to  be  determined  to  ensure  that  the  integral  mass  flux  relation 

/  pU’^dA  =  Q  =  constant  (3.10) 

A 


is  satisfied.  In  the  present  work,  pm  is  adjusted  by  a  secant  iteration 
during  solution  of  the  streamwise  momentum  equation  (2.8)  so  that  integral  mass 
flux  relation  (3.10)  is  satisfied. 

3. 3  Summa  rv 


A  summary  of  the  procedure  used  to  advance  the  solution  a  single 
streamwise  step  to  the  (n+1)  level  from  previously  computed  quantities  at  the 
n  level  follows.  Coordinate  and  vector  basis  information  required  to  solve  the 
primary /secondary  equations  are  prescribed  input.  Governing  equations  are 
assumed  to  be  linearized  and  discretized  as  previously  described. 

(1)  Equations  (2.9)  and  (2.10)  form  a  linear  coupled  system  for  i|in+l 

and  which  is  solved  as  a  (2x2)  coupled  system.  For  this 

purpose,  artificial  time  derivatives  are  added  to  each  equation  and 
an  iterative  block-implicit  scheme  [2]  is  used.  In  prescribing 
no-slip  boundary  conditions,  the  tangential  component  (3.8)  contains 
a  contribution  from  ;  this  contribution  is  evaluated  using  <f> u . 
Linearization  (3.4)  provides  that  terms  involving  other  unknown 
dependent  variables  are  evaluated  at  the  n-level. 

(2)  The  pressure  equation  (2.11)  is  solved  for  using  an  iterative 

scalar  ADI  scheme.  Updated  values  of  vJJ+1  and  are  available 

while  other  dependent  variables  are  evaluated  at  the  n-level. 

(3)  The  energy  equation  (2.12)  is  solved  for  En+*  using  a  scalar  ADI 
scheme . 
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(4a)  Solution  of  the  streamwise  momentum  equation  is  described  for 

internal  flow.  Using  an  assumed  value  of  pn+*  to  begin  a  secant 

m 

iteration,  the  streamwise  momentum  equation  (2.8)  is  solved  for 
u^+1  using  a  scalar  ADI  scheme.  The  streamwise  pressure  gradient 
is  evaluated  at  the  (n+1)  level. 

(4b)  The  density  pn+l  is  evaluated  algebraically  from  the  equation  of 
state  (2.13) 

(4c)  The  integral  mass  flux  through  the  transverse  plane  is  computed 
using  (3.10). 

(4d)  Assuming  the  initial  guess  for  was  not  exact,  the  integral 

mass  flux  relation  (3.10)  will  not  be  satisfied.  Steps  (4a-c)  are 
repeated  iteratively  using  a  standard  secant  method  [3]  to  find  the 
value  of  which  leads  to  u^+l  and  ph+l  satisfying  the 

integral  mass  flux  relation  (3.10).  In  the  incompressible  limit, 
all  equations  in  the  secant  iteration  are  linear  and  the  third 
iteration  is  exact.  At  high  subsonic  Mach  numbers,  a  fourth  or 
fifth  iteration  may  be  required. 

(5)  Finally,  the  scalar  potential  equation  (2.7)  is  solved  for  i})n+l 
using  an  iterative  scalar  ADI  scheme. 
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MOTION  OF  SHROUD 


TIP  CLEARANCE  GAP  =  0.01  CHORD 

INBOARD  AND  LATERAL  BOUNDARIES  ARE  0.5  CHORDS  FROM  TIP 
Re  =  20,000  (LAMINAR  FLOW) 

SHROUD  VELOCITY  =  0.5  PRIMARY  FLOW  VELOCITY 


Figure  A.  Schematic  of  Cross-Section  Normal  to  Primary  Flow 
Direction  for  Shrouded  Tip  Flow  Application. 
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Figure  la.  Single  Airfoil. 
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Figure  lb.  Zero  Pitch  Planar  Cascade. 


Sheared  Airfoil  Coordinates  Blade,  Thickness  from  Eq.  1 
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